Summary of Research 


NAG3-2857 


Improved Temperature Dynamic Model of Turbine Subcomponents 
for Facilitation of Generalized Tip Clearance Control 

Javier A. Kypuros, Ph.D., Rodrigo Colson, and Affedo Munoz 
University of Texas-Pan American, Edinburg, Texas 


Abstract 

This paper describes efforts conducted to improve dynamic temperature estimations of a turbine tip 
clearance system to facilitate design of a generalized tip clearance controller. This work builds upon 
research previously conducted and presented in [1] and focuses primarily on improving dynamic 
temperature estimations of the primary components affecting tip clearance (i.e. the rotor, blades, and 
casing/shroud). The temperature profiles estimated by the previous model iteration, specifically for the 
rotor and blades, were found to be inaccurate and, more importantly, insufficient to facilitate controller 
design. Some assumptions made to facilitate the previous results were not valid, and thus improvements 
are presented here to better match the physical reality. As will be shown, the improved temperature sub- 
models, match a commercially validated model and are sufficiently simplified to aid in controller design. 

Nomenclature 


Symbol 

Units 

Description 

A 

_2 

m 

surface area 

V 

3 

m 

volume 

L 

m 

characteristic length (L = V / A) 

T 

°C 

temperature 

e 

°C 

temperature differential {6 = AT = T(x,() — T(x,t = 0) = T -T) 

c 

J/kg-°C 

specific heat 

h 

W/m 2 -°C 

convection heat transfer coefficient 

k 

W/m-°C 

thermal conductivity 

r 

m 

radius 

t 

sec 

time 

w 

m 

width 

a 

m 2 /sec 

thermal diffusivity (a - k 1 pc) 

S 


penetration depth of surface thermal effects 

P 

kg/m 3 

density 

V 


Poisson’s ratio 

V 


film cooling effectiveness 


1 Introduction 

The focus of this study is to improve estimations of the temperature distributions in the rotor, blades, 
and casing/shroud substructure. Results from the previous model were compared to a validated 
commercial model. Temperature estimations for the rotor and blades were found to be erroneous, and the 
casing/shroud substructure was modeled differently. Since the intent of this model is to sufficiently 
estimate dynamic turbine tip clearance to facilitate controller design, temperature estimations need to be 
improved while maintaining a relatively simple model. The major culprits responsible for errors are 
suspected to be a few of the assumptions used to assist in deriving the previous analytical solutions. 

The commonly-used semi-infinite solution to the conduction equation with convective boundary 
condition was previously employed to estimate the surface temperature of the rotor, blades, and shroud: 
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pck 
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where T w =T(x = 0,1) , T, =T(x,t = 0) , and T x are the surface, initial, and fluid temperatures. This 
solution was derived making the following assumptions: 

1 . the heat flows primarily one-dimensionally, 

2. the surface is exposed to convection, 

3. the temperature of the convective fluid is constant, 

4. the average temperature of the component is approximately equal to the surface temperature, and 

5 . the material properties are assumed to be constant. 

The first two assumptions are reasonable but the remaining three are questionable. The temperatures of 
the fluids flowing over respective surfaces vary with time and that must be accounted for. Some of the 
components are sufficiently massive that their thermal inertia must be considered. The temperatures vary 
over significant ranges suggesting that the material property temperature-dependencies may need to be 
considered. These assumptions are critically evaluated to improve temperature estimations. 

Additionally, the steady-state conduction equation applied in the radial dimension was used to 
determine the temperature distribution through the shroud so that the thermal stress could be calculated, 

T(r) = T + (T -T ) ^L/Jjnnsrl n\ 

V / x S, inner ' V s,outer * s, inner J i / / \ V~/ 

^y outer ' ^ inner ) 

where T Siinner and T St0Uter are the inner and outer surface temperatures, respectively. To facilitate this 
solution, the heat transfer through the shroud was assumed to be quasi-steady-state (i.e. the heat transfer 
in the shroud was assumed to occur relatively rapidly). As mentioned above, the inner and outer surface 
temperatures were calculated using the semi-infinite solution. Thus, the same questionable assumptions 
apply. Furthermore, the quasi-steady-state assumption was never fully validated. 

In the subsequent sections, alternative formulations are employed to improve estimation of the 
temperature distributions. The three major methods used to gamer analytical solutions to the transient 
conduction equation are the lumped-capacitance method, the Laplace transform method, and the 
approximate method of assumed solutions. The remainder of this document describes the resulting 
solutions and shows which methods provided the best estimations for each of the components. 


2 Improved Temperature Estimation 

The methods employed here to derive analytical solutions are intended to critically evaluate the 
validity of the questionable assumptions detailed above. The rotor is first evaluated and alternative 
solutions that incorporate ransient boundary conditions are derived. As will be shown, the transient 
boundary conditions are more important in approximating the average temperature of the rotor than the 
temperature-dependent material properties. Further investigation of the turbine blades suggested that the 
relatively low ratio of volume to surface area of the blades would facilitate use of the lumped-capacity 
method. Results detailed below will show that the average blade temperature can be reasonably 
approximated by assuming constant material properties and incorporating the transient boundary 
condition into the lumped-capacity method. The casing/shroud substructure is the least investigated in the 
literature so a formulation similar to that suggested for the rotor is detailed below. However, at this time 
there are no readily available validated results to confirm the accuracy or appraise estimation 
improvement due to this formulation. 
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2.1 Rotor 

As previously mentioned, the rotor temperature was originally calculated using Equation (1). 
Equation ( 1 ) was derived using the Laplace transform method and assumes that the convective surface is 
exposed to a fluid at a constant temperature. Additionally, the average temperature was assumed to be 
approximately equal to the surface temperature. To better calculate average rotor temperatures, the 
temperature distribution should be integrated with respect to x (the distance measured from the surface 
across the width of the solid) and divided by the rotor width, w, 


e. 


average 


i l 

= — f 0(x,t)dx. 

W f 


By applying Laplace transforms 1 , to the 1-D conduction equation, 

39 3 2 9 

— = oc — r, 

dt dx 2 


( 3 ) 


( 4 ) 


it can be shown that the temperature distribution in a semi-infinite solid is 
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where 9 m =T a> -T r Note that for a constant fluid temperature (i.e. 6 m - constant) the above solution 

evaluated at the surface (x = 0) simplifies to Equation (1). The presence of a convolution (*) and the 
multiplication of an exponential and error function make this equation difficult to integrate analytically 
with respect to x to calculate the average temperature. Moreover, as shown in Figure 1, the rotor does 
not have only a single surface exposed to convection as is assumed for a semi-infinite solid but two 
surfaces. 


insulated 

(assumption) ^ 



Figure 1. Rotor schematic: (a) frontal view and (b) cross-sectional view. 

To facilitate a simpler, closed-form approximation that can accommodate the two convection 
boundary conditions, the method of assumed solutions is employed. This approximate method as detailed 
in [2] utilizes the integral form of the conduction equation and an assumed temperature profile (usually 
defined as a polynomial which is easily integrated) that matches the boundary and initial conditions. The 
integral form of the conduction equation is [4,5] 


1 Refer to [2,3] for use of the Laplace transform in solving the 1-D conduction equation. These references also 

include useful transform tables. 
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where 8 is the penetration depth shown in 


Figure 2. It is the depth from the surface to which the thermal effects are felt (i.e. the distance from the 
surface where the temperature gradient goes to zero). 


i 



T(w,t) 




Figure 2. (a) Infinite slab temperature distribution and (b) insulated centerline assumption. 

As shown in Figure 3, as t — » oo the thermal effect of the convection fully penetrates (i.e. S = w/2 ) 
at which point the temperature profile is fully developed and begins to move up the temperature axis. 
After this, the penetration depth remains constant. Prior to this instant the temperature change at the 
centerline remains zero. Hence, to properly capture the temperature profile, the solution must be broken 
into two parts: (1) for S <wf 2 and (2) for 8 - w/2 . Moreover, while 8 <w/ 2 , the temperature 
profile only applies to the interval 0 < x <8 because over the interval 8 < x < w/2 the temperature 
differential is zero (i.e. 8(x,t) = 0 for 8 < x < w/2 while 8 < w/2). 



Figure 3. Dynamic temperature distribution in an infinite slab. 

To approximate the average temperature of the rotor, the temperature distribution through the width 
of the rotor must be determined. The distribution at a radius sufficiently far away from the outer 
circumference can be approximated by treating the cross section of the rotor like an infinite slab with 
convection at both surfaces. It is assumed that both surfaces are exposed to the same fluid temperature. 
Thus the profile through the width will be symmetric about the centerline as shown in Figures 

Figure 2 and Figure 3. Since the distribution is symmetric the half-slab case can be used (refer to 
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Figure 2 (b)) and the temperature can be averaged over half the width of the rotor. The boundary 
conditions then are, 


at x = 0 


at x = w/2 


(±t, A 

K'~ 


dx 


= 0 


'±e. 


k 

k, 

\ 




and 


dx 


= 0 . 


y x=v! 2 


Though the boundary conditions can be satisfied by using a quadratic distribution, significantly more 
accurate results are garnered by using the following cubic distribution [2], 


-0 compressor (t)[\ - 


Sit), 


(7) 


By evaluating the above equation at x = 0 one can solve for 0(0 J) , 


q~(0.0= -- y.. ft 


hS + 3 k r 


compressor * 


Substituting this into Equation (7) gives 


S rotor i^’ t) 


K5 


- 0 . .. 


h r S + 3k r compressor 

which has a partial derivative with respect to x of. 


1 -- 
V S j 


( 8 ) 


L e rLL e 

W rotor }j S ~\~3k U(:om P ressor 


f 

1 -- 
V S J 


(9) 


The integral form of the conduction equation (Equation (6)) is used to calculate 8 . Due to the 
boundary conditions, the second term on the left side and first term on the right of that equation are zero 
(i.e. 0 rolor (8,t ) = 0 and {dO rolor / dx) x=s = 0 ). Thus using the assumed temperature profile in Equation 
(8) and the temperature gradient in Equation (9), Equation (6) simplifies to 


d_ 

dt 


KS 

h r 8 + 3 k r 


0, 


compressor 


SJ 

it 1 


Y 

~ dx 

3 h r a 

dx 

J 

~ ^rotoriS *t)~~ ~ OC f 

dt 

_h r 8 + 3k r 


which can be solved for 8 , 


8 = 


\2a r {h r 8 +3k r ) 


ll/2 


compressor 


6 }h8+3k 


dr 


compressor 


0 r^ 


( 10 ) 


Equations (7) and (10) can be used to simulate the dynamic temperature distribution up until the 
penetration depth reaches half the width of the rotor (i.e. while 8 <w/2). Once the penetration depth 
reaches half the width of the rotor (i.e. 8 = w/2), the thermal effects of the surface convection have 
penetrated fully and the temperature profile begins to rise. The half-slab temperature profile continues to 
be approximately cubic but now the distribution is 
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0ro<or(X,t) = [0 r coX^)-9 wn (t)i\ +O w/2 {t) 

\ w 1 2 ) 


(ii) 


where @ w/2 (t) is the centerline temperature differential which is zero until the penetration depth reaches 

the centerline. By evaluating the gradient of the assumed temperature profile at the surface and applying 
the convection boundary condition it can be shown that 


® rotor (0? 0 = 


h r (w/2)0 CO m pressO r+3k r 0 wn 


h r (w/ 2) + 3k r 

which when applied to Equation (11) gives the following temperature profile 


e rotoA x ’ t ) = 


hfw/2) 
h r ( w/2) + 3k r 
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w/2 


+e. 


w/2 


( 12 ) 


and thermal gradient 


dx 
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By using Equations (6), (12), and (13) the rate change of centerline temperature can be calculated as 


@w!2 ~ 
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2>h r {w!2) + \2k r 


-h r (w/ 2)B 
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compressor 
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(13) 


(14) 


which can be integrated to determine the centerline temperature. Equations (12) and (14) can be used to 
simulate the half-slab temperature profile after the penetration depth reaches the centerline. 


To summarize, the general equation for the approximate temperature distribution is 


0 rotor CM) 


hS 


hS + 3k. 


compressor ^w/2 ) ^ 
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(15) 


where 0 w{2 = 0 for S <w/2. The penetration depth can be calculated using Equation (10) until it 

reaches the centerline at which point it remains constant and equal to half the rotor width. Thereafter, the 
centerline temperature is calculated using Equation (14). Equation (15) can be substituted into Equation 
(3) to calculate the width -wise average rotor temperature, 


(O, 


+{3h r s -nk.y>. 


w/2 


average 


4h r S + \2k r 


(16) 


The rotor width will likely vary over the radius so the average width can be used to provide a reasonable 
approximation using the above equations. 


2.2 Blades 

As has already been mentioned for the turbine blades, the low ratio of volume to surface area suggests 
that the lumped-capacity method is applicable. More specifically, for most turbine blades the Biot 
number is small, Bi = hL / k « 1 . Therefore, the approximate average blade temperature is 
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@ blade ~ blade 0 ) 


(17) 


where L is the characteristic length ( L — V / A ) and the reference temperature differential, 6 ref , is 
calculated using the film cooling effectiveness. 


77 = 


^ref ^ turbine 

a j a 

^compressor ^ turbine 


(18) 


(Note V and A are the blade volume and surface area.) As shown in Figure 4 the mainstream flow and 
film cooling flows have temperatures T lurbine and T compressor , respectively. Equations (17) and (18) are 
used to simulate the dynamic blade temperature. 



film cooling hole 


3 


T blade & h 


blade surface 


(b) 


. Figure 4. Blade schematic: (a) turbine blade and (b) blade surface. 


2.3 Shroud/Casing 

The shroud is a series of arcs circumferentially surrounding the turbine rotor and blades (refer to 
Figure 5). The arcs are mounted to the casing. Since the shroud is not one continuous rigid structure it is 
believed that the casing to which the shroud is attached dictates the radial deflection of the shroud. The 
inner surface of the casing is exposed to compressor discharge air. The outer surface conditions are 
questionable. To facilitate an approximation, it is assumed that turbine is mounted in a test rig where the 
outer casing surface is exposed to the atmosphere. 
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One could approximate the dynamic temperature distribution in the casing by using the cylindrical 
form of the conduction equation. However, this may unnecessarily complicate the solution because the 
wall thickness of the casing may be several orders of magnitude less than the radius of the casing (i.e. 
(r 2 -r,)/r 2 « 1 ). Therefore, the casing can be approximated as a plane wall with two distinct boundary 
conditions as shown in Figure 6. This simplifies analysis and enables the use of a formulation similar to 
that used in §2.1. By using the principle of superposition, the two-boundary-condition problem can be 
broken in to two problems: (1) the temperature distribution due to thermal effects at the inner surface of 
the casing and (2) the temperature distribution due to the thermal effects at the outer surface of the casing 
(refer to Figure 6). One can calculate the temperature distributions due to convection at each surface and 
then sum the two solutions to approximate the overall distribution. Once the thermal effects of each side 
penetrate the width, the fully developed profile moves up the temperature axis as in Figure 7. 


7 aim P at 


& Pam, 




Figure 6. Infinite slab with two distinct boundary conditions: (a) temperature distribution, (b) distribution due to 

7 1>q0 and (c) distribution due to r 2)0 ©. 
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Figure 7. Temperature distribution in an infinite slab with two distinct convection boundary conditions. 

Since the thickness of the casing is significantly less than the width of the rotor it is believed that a 
quadratic distribution (as opposed to a cubic distribution like that used for the rotor) should provide a 
sufficient approximation of the temperature distribution. The casing solution must be defined for the 
interval that the penetration depth is less that the casing thickness (i.e. 8 < w) and for the interval when it 
is equal to the casing thickness (i.e. 8 = w ). To simplify derivation of each distribution, two coordinates 
are used: x\ and x 2 . The first is measured from the inner surface edge and second from the outer surface 
edge. The two distributions can be easily summed by recognizing that x 2 = w - x, . The temperature 
distributions can be derived in a similar fashion to the way the rotor distribution was. 

By using a quadratic distribution and applying boundary conditions similar to those in Figure 6(b), 
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it can be shown that the approximate temperature distribution due to the compressor discharge air is 
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compressor 






(19) 


where 6 X v , the outer surface temperature differential due to convection at the inner surface, is zero until 
the inner penetration depth reaches the outer surface (i.e. 8 X = w ) at which point it is calculated using 




1 

2h l w+6k c 


— hxvf) — f) ) 

rl\ WU c(jm pr e ssQr • ^ \y compressor C7 1 ,w / 


( 20 ) 


which is derived using Equations (6) and (19). One can also use Equations (6) and (19) to derive the 
penetration depth due to thermal effects at the inner surface prior to reaching the outer surface, 


< y ,= 




-, 1/2 
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compressor 
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dr 
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The component of the average temperature differential due to 0 ] is 

fo), 


compressor + {^ ] S ] + 6k c )# lw 


/ average 


3f\S } +6k c 


(22) 


Similarly, one can determine the temperature distribution due to convection to atmosphere at the outer 
surface in terms of coordinate x 2 . Assuming for simplicity that the atmosphere is at a constant 
temperature, and that the boundary conditions are like those in Figure 6(c), 


at x 2 = 0 


at x 2 = w 


0 2 (S 2 ,t) = 0 


{dx ) 

( d > 


x 7 =0 
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dx 
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the temperature distribution due to the thermal effects at the outer surface is 
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0 2 (x 2 ,t) = 


, „ (0--O 

h 2 w + 2k c \ w j 


0 . 


2,w 


(23) 


where 0 2 w , the inner surface temperature differential due to convection at the outer surface, is zero until 
the outer penetration depth reaches the inner surface (i.e. S 2 = w), and then it is calculated from 


a 3a e /z 2 

2,w w(h 2 w + 3k c ) 


(@atm @2,w ) 


(24) 


which is also derived using the integral form of the conduction equation. The component of the average 
temperature due to 0 2 is 


te), 


MA. + (ih,o 2 + 6k c Y> 


2,v 


average 


3h 2 S 2 +6 k c 


(25) 


The overall average temperature can be approximated by summing Equations (22) and (25). 


3 Results 

Based on results for a validated model, it was determined that the previous turbine tip clearance 
model erroneously estimated the temperatures for the shroud/casing, blades, and rotor. As shown in 
Figure 8 the commercial model ( — CES) estimated a slower transient for the rotor and a more pronounced 
transient for the blade relative to the original turbine tip clearance model ( — ACC). 
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Figure 8. Comparison of original model (ACC) to validated results (CES). 

By implementing the corrections detailed in §2. land §2.2 the new simulation produces the results in 
Figure 9. Note how the new results ( — ACC2) for the rotor almost exactly match the CES estimation. 
The blade results are also much improved over the previous model. Furthermore, as was detailed in [1] 
many of the parameters were estimated based on numbers found in literature, and with minor adjustments 
or improvements to the values used, the new model can match the CES model more closely. 



T 1 1 1 1 1 * 1 



Figure 9. Temperature estimations using approximations detailed in §2.1 and §2.2. 

Figure 1 0 shows the average temperature estimation for the casing using the analysis from §2.3. The 
results are similar in form to the estimation for the rotor. 
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Figure 10. Temperature estimation for casing using approximation in §2.3. 


4 Discussion 

As was expected, some of the assumptions used in the original model were inappropriate, particularly 
assumptions 3 and 4. Assumption 5 appears to be reasonable because the improved results detailed above 
were achieved without incorporating temperature-dependent material properties. However, one should 
note that the approximate method used for the rotor and casing can be adapted to incorporate material 
property temperature dependencies. By using a transformation suggested by Goodman [4], 

e 

v(x,t)= J/x dd , 

0 


in the conduction equation (Equation (4)), one will yield 
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which can be integrated with respect to x for an equivalent integral form of the conduction equation, 
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By using the approximate method of assumed solutions, a estimation incorporating temperature 
dependencies could be derived with some work, but based on the current results the added effort would 
likely be unnecessary. 

In conclusion, the improved estimations provided herein would substantially advance the existing 
model while maintaining the relative simplicity desired to facilitate design of the future active clearance 
control. 
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